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Abstract 

We present the SIMUB package developed at Dubna for MC 
generation of B meson production and decays. The starting version 
of the package includes lepton modes of B decays, in particular, 
semileptonic decays B D^*Hv and "golden" mode B J/ijj{-^ 
/i/z) KK) with taking into account all theoretical refinements 
including B^ — B^ oscillations and angular correlations. 

SIMUB is a Monte Carlo (MC) generator of 5-meson production and 
decays which is under developing at Dubna for the Compact Muon Solenoid 
(CMS) Project at CERN (see SIMUB documentation in [1]). The main 
motivation for this activity was that already existing generators do not 
take into account the theoretical refinements which are of great impor- 
tance for MC studies of -B-decay dynamics. In particular, in the gen- 
erators PYTHIA [2], QQ [3], and EvtGen [4], the time- dependent spin- 
angular correlations between the final-state particles are not included in 
the proper way for the so called "golden" decay B^(t), Bg{t) J/'^{-^ 

H~) (/){—>■ K^K~). The dynamics of this decay is described by four- 
dimension probability distribution function depending on decay time and 
three physical angles. The algorithms of multidimensional random num- 
ber generation have been elaborated and then implemented in the package 
SIMUB to provide tools for MC simulation of sequential two-body decays 
B^{t) ,B{t)^ 0102) &(— * &1&2) in accordance with theoretical time- 
dependent angular distributions. 

This paper is organized as follows. First, some general theoretical as- 
pects are considered in the context of B^ — B oscilations and angular 
correlations. Then, we discuss the algorithms of multidimensional ran- 
dom number generation restricting ourself by detail consideration of only 
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B Vi{^ n+n-)V2{-^ P+P-) channel and its J/ip(p "golden" 

mode. Finally, we present a general information about SIMUB package 
and compare time and angular distributions for "golden" decay mode gen- 
erated by SIMUB with those obtained by using PYTHIA. 

1. B'^-W mixing and time evolution of neutral 
5-meson states 

The time dependence oi B^ ^ f decays, in which / is not a CP eigen- 
state, is not purely exponent due to the presence oi B^ — B mixing. This 
mixing arises due to either mass difference or decay-width difference be- 
tween the mass eigenstates of the B^ — B system. The time evolution of 
the state of initially, i.e. at time t = 0, present B^ (^) 

meson can be described as follows: 



B'{t)) ^ g4t)\B') + g4t)\-^) , = 0) = 1 , ^_(i = 0)=0; 

^{t))^g4t)\B')+g_{t)\'^), g^t ^0)^0, g_{t ^ 0) ^ 1 . 



The strong eigenstates \B^{t)) and \B (t)) are not eigenstates of the 
full Hamiltonian due to the weak interaction. Thus, they most generally 
evolve according to the Schroedinger equation: 



Diagonalization of full Hamiltonian M — (i/2)r (see [5] for more detail) 
gives 



Here fi± = Mljh — {i/2)VLiH are eigenvalues of the full Hamiltonian cor- 
responding to the masses and total widths of physical ''lighf and "/ieav?/' 
eigenstates of full Hamiltonian \Bl/h) = \B±), and 






g+it)=g-{t)/a\ 



g.{t)=g+{t). 




(1) 
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is a phase factor defining CP transformation of flavor eigenstates of neutral 
5-meson system: CP\B^{t)) = r]^p\^{t)). 

In Eq. (1) we have made a good approximation <S M12 valid for 
S° — 5° system. This approximation implies that there is no CP violation 
m. — B mixing, i.e. the probability for B^ to oscillate to a S is equal 
to the probability of a to oscillate to a 

Prob(S° ^ S°) ^ Prob(S° ^ 

Such an asymmetry in mixing, which results from | a | 7^ 1 , is often reffered 
to as indirect CP violation. The situation with indirect CP-violation in 
B^ — B mixing is in contrast to the — K system, in which 



Prob(ir ^ K"") > Prob(i^'° ^ K^). 

The indirect CP- violation is the main source of CP- violation in neutral 

K decays. 

From CPT theorem it follows that 

Mn = M22 = M , M12 = M*i ; 

Pii = r22 = r , ri2 = r2i , 

where M and F are "mean" values: M = {Ml + Mh)/2, F = (Fi + FH)/2. 
The values AM = (Mh - Ml) and AF = (Fj:^ - Tl) are used also. 

The explicit expressions for mixing probabilities as a function of time 
are given by 

Prob(P° at t\B° at i = 0) = Prob(P° at t\B^ at i = 0) 

^ [e-^^* + e-^«* + 2e-^*cos(AMi)l , 
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Prob(P° at t\B° at i = 0) = ^|a|2[e-^^* + e'^^* - 2e-^*cos(AMi)] , 
Prob(P° at t\B" at i = 0) = [e"^^* + e"^"* - 2e-^*cos(AMi)] , 

(2) 

For case |q;| = 1 (no indirect CP violation) in approximation F^ ~ Fh, 
the latter two expressions in Eqs. (2) can be rewritten in form 

Prob(# at at t = 0) = Prob(5° at t\B^ at t = 0) ^ sin^ 



which coincides with the formula for — B mixing probabihty used in 
PYTHIA [2]: 

Vfup = sin^ j . (3) 

Here (r) = l/F is the mean hfetime, and x is the mixing parameter [6]: 

x = Xd = AM/r = 0.74 ± 0.02 in 5^ - B^ system; 

x = Xs = AM/r > 19 in 5° - B^ system. 

In PYTHIA, the initial B^ meson is allowed with the probability (3) to 
decay like a and vice versa. 

2. Time evolution of decay amplitudes 

Time evolution of the amplitudes of transitions \B^[t)),^{t)) — >■ |/) 
induced by Hamiltonian H^ff is represented by 

Af{t) = {f\neff\B\t)) = g+it){f\neff\B')+g4t){f\neff\B') , 
Af{t) ^ = g+m\neffm+g4t){f\neff\'^) . (4) 

If I/) is the eigenstate of CP-operator 

CP\f) = r?^p|/) , {4p = ±1) , 
the relations (4) can be rewritten in form 

L VcpVcp • 

Af{t) = iflHefflB') \-g+{t) + 9-{t)^^Cf 
L VcpVcp ■ 

Here 

^' ~ {miT,m • 

where Ti-^ff = CPTCeffCP, and 0/ is the CP-violating weak phase [7]. 

The CP- violating weak phase (pf is introduced through interference ef- 
fects between B^ — B^ mixing and decay process with final state |/) being 
the eigenstate of CP operator. In case of decays 5°, B^ —>■ f, where 



q & d,s, the value of weak phase ^^^^ can be expressed in terms of ma- 
trix elements of the combinations Qp'' of four-quark operators and Wil- 
son coefficients involved in the low energy effective Hamiltonian i^e// (see 
Refs. [7, 8]): 



(5) 



Here 0^^^ is the — mixing phase: 

where Vij is the Cabibbo-Kobayashi-Maskawa (CKM) matrix elements 
while P and ^7 are standard parameters related with angles of the uni- 
tary triangles. 

In general, the observable ^j:*^ suffers from large theoretical uncertain- 
ties in the hadronic matrix elements {f\Q''^\'B^). However, if the decays 
B^,^g ^ f are dominated by a single CKM amplitude, the corresponding 
matrix elements in Eq. (5) are canceled, and ^j?^ takes the simple form: 

where q & d, s, j & u, c, and 0;^g^ is a CP-violating weak decay phase: 

, J _ /y* Y \ _ \ for j = u - dominant b — > uuq , 

''^'^ ~ •'^ ■'^ 1 for j = c ~ dominant b — > ccq , 

where 7 is the angle of unitary triangle. 

3. Transversity amplitudes and time-dependent 
observables for decays B{t) — ^ V1V2 

The decay B^{t) — > V1V2, where Vi and V2 are vector mesons in final 
state, is described in terms of three transverse amplitudes Aq, A\\ and A± 
[9]: 

(6) 
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Here Ey^ is the energy of the V2 in the Vi rest frame; p is a unit vector in 
the direction of the momentum of V2 in Vi rest frame; eVi and eVa are the 
polarization three- vectors in the V2 rest frame; = p - e is the longitudinal 
component and is the transverse component of the polarization vector. 
The decay B (t) V1V2 is described in analogous way in terms of the 
transvcrsity amplitudes Aq, A\\ and Aj_. 

Eq. (6) corresponds to decomposition of final state of — > V1V2 de- 
cay over transversity basis which can be represented in form of 3-dimensional 
vector 

|/) = (|(Vl^)o),|(V^1^2)||),|(W2U)). 

The transversity components are eigenstates of CP-operator: 

CP\{V,V2)f) = vtpmv2)f) , (/ = 0, II, ±) 
with eigenvalues 

Vcp — 1 ' Vcp — 1 ' Vcp — ~1 ■ 

Thus, the Eq. (5) with (/|7ie//|-B°) = Af{0) can be applied to calculate 
the time evolution of transversity amplitudes Af{t) (/ = 0, ||, _L). 

Matrix clement squared involves into consideration the following six 
time-dependent bilinear combinations of transversity amplitudes {observ- 
ables): 

0,{t) = |Ao(t)p . 02{t) = , 0,{t) = \A^{t)\^ , 

0,{t) = lm(Al{t)A^{t)) , 0,{t) = Re{Al{t)A\\{t)) , 

0,{t) = hn{Al{t)A^{t)) , (7) 

in case of decays B^{t) V1V2, and similar combinations of Af{t) (/ — 
0, II, ±) for decays B\t) V1V2. 

4. The time-dependent angular distributions 

Let us consider sequential two-body decays 

, ^ aias) 6162) , 

where a and b are vector mesons decaying into pairs of particles ai, 02 and 
61, 62, respectively. The angular distributions for these decays are governed 
by spin-angular correlations (see [10]- [12]) and involve three physical an- 
gles. In case of the so-called helicity frame [11] these angles are defined as 
shown in Fig. (1). The 2;-axis is defined to be the direction of 6-particle in 
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the rest frame of the B'^. The x-axis is defined as any arbitrary fixed di- 
rection in the plane normal to the z-axis. The y-axis is then fixed uniquely 
via y — X X z. The angles (©oi, XaJ specify the direction of the oi in the 
a rest frame while {Qbi, Xbi) ^-re the direction of bi in the b rest frame. 
Since the orientation of x-axis is a matter of convention, only the difference 
X — Xai ~ Xbi of the two azimuthal angles is physical. 




Figure 1: The decay B — > a(— > 0102) b{-^ bib2) in helicity frame 

In most general form the angular distribution for decay B^{t) 
0102) > &162) can be expressed as 

d^T 9 ^ 

(8) 

where Oi {i — 1, 6) are the observables (7) and gi are the functions of 
physical angles Qai, ©61 and x- 

The time- dependent angular distribution function (8) is used as density 
of the probability function for MC simulation of the vertex and kinemat- 
ics of the final-state particles in case of decay -B°(t) — > 0102) 
6162)- The variables of the function /o(©ai, ©61, X; ^) can not be fac- 
torized and randomly generated in the independent way. Nevertheless, 
the random generation of cos 6^ , cos 6^^ , x and t can be performed ei- 
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ther simultaneously according to distribution function (8) by using four- 
dimensional random generator or successively, one after another, by using 
single-dimensional random number generators with accordance to distri- 
bution functions obtained by successive integration of the function (8) over 
its variables. 

Let as consider the latter approach in the case of sequential random 
generation of the variables t, x, cosQai and cosO^^. The following three 
distribution function are used in this case: 

r+l 

/i(0ai,x;t) = dcosQbiMQai,Qb„x;t) , 

/2(x; t) = dcosOa, fi{Qa„x; t) , (9) 

and 

Mt) ^ rdxf2{X-:t) = 0,{t) + O^it) + O^it) . (10) 

Jo 

The procedure is as follows: 

• First, the proper time t is randomly generated according to the dis- 
tribution function f-i{t). 

• Second, the angle % is randomly generated according to the single- 
dimensional distribution given by function /2(x; t) with t being fixed 
to be equal to the time, generated at the first step. 

• Then, the value cos Qai is generated according to distribTition func- 
tion /i(Oai,X;^) with values of t and x fixed to be equal to their 
values generated at the previous two steps. 

• Finally, the value of cos 0^^ can be generated according to the single- 
dimensional distribution given by function /o(0an 06i) Xi ^) with a 
properly fixed values of x and cos0ai. 

The functions of Eqs. (9) and (10) present the most important experimen- 
tally observable distributions. 

5. Decays B^{t),B^{t) Vi{l+l-)V2{P+P^) 

In case of sequential decays B°{t) ,^{t) ^ Vi{^ /+/-)\/2(^ P^P~), 
where and are lepton and pseudoscalar mesons, respectively, let 
associate the particles in the final states with general notations used in 
Fig. 1 as follows: 

Vi ^ a , Oi , /~ — >■ 02 , 

V2^b, P^h, P'^b2. 
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In this case the functions gi{Qai, 0fei,x) Eq. (8) are defined as [11]: 
gi = 4sin^0aiCos^06j , 

512 = (1 + cos^OaJsin^Ob, - sin^eaiSin^e6iCOs2x, 
grg = (1 + cos^eai)sin^©6i + sin^©„iSin^e6iCOs2x , 
g4 = 2sin^6a,sin^O{,iSin2x , 
t/5 = -V2sin20a^sin206jcosx, 

ge = V^sin2©a^sin2©;,jsinx . (11) 
The exphcit form of the distribution functions (9) are given below: 
/i(cose„,,X;t) = ^Oi(t)sin2e„, + ^02(t)(l-sin2e„,cos2x) 

+^03{t) (1 - sin'e„,sin\) - O^it) sin'e„iSin2x , 
Mx,t) = 0,{t) + ^02{t) {3 -2cos'x) 

+ ^03{t) (3 - 2sin2x) - 04t) sin2x . (12) 

6. Physics parameters for decay J/ '^4' 

For numerical MC simulation of the decay J/'^4' the following 

physics parameters should be defined: 

• q; = 1 - no CP violation in — 5^-mixing; 

• Ml = 5.3696 GeV - mass of the "light" B meson; 

• AM = Mh-Ml = 10.6 ps-i; 

• V I = \/t - width of the "light" B meson, where r = 1.493 ps; 
. AT = Th-Tl = 0.21^; 

• (t)f^ — 0.03 -i- 0.06 - CP- violating weak phase; 

• initial values of observables at time t = 0: |Ao(0)|, |A||(0)|, |y4_L(0)|, 
and two CP-conserving strong phases 5i = arg[74|j(0)74j_(0)] and 82 = 
argK(0)A^(0)]. 

The time-reversal invariance and naive factorization lead to the follow- 
ing common property: 

Im[^S(0)^±(0)] =0, lm[^f|(0)^x(0)] - 0, 
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ReK(0)A||(0)] = ±|Ao(0)A||(0)|. 

Moreover, in absence of strong final-state interactions, 5i = vr and 82 = 0. 

Using low-energy efi^ective Hamiltonian [8] for nonlcptonic 6-quark tran- 
sitions b — > sec, the amplitudes Af{0) (/ = 0, ||, ±) of J/i'4> decay 
can be calculated within the framework of naive factorization in terms of 
the form factors of transitions B ^ (j) induced by quark currents. The 
B ^ (f) form factors can be related to the B K* case by using SU(3) 
flavor symmetry. In table 1 we collect the predictions of Ref. [13] for the 
ratios of S° — > J/ip(t) observables calculated with B ^ K* form factors 
given by different models [14]-[16]. The quantity 



ro(o) 




ro(o) + rT(o) 


Ao(0)|2 + |A)|(0)|2 + 





describes the ratio of the longitudinal to the total rate at t = 0. 



Table 1. Predictions for B^ —>■ J/ip4> observables 



Observable 


BSW [14] 


Soares [15] 


Cheng [16] 


l^||(0)|/|Ao(0)| 


0.81 (0.77) 


0.82 (0.78) 


0.75 (0.70) 


l^±(0)|/|^(0)| 


0.41 (0.40) 


0.89 (0.88) 


0.55 (0.54) 


ro(o)/(ro(o) + rr(o)) 


0.55 (0.57) 


0.40 (0.42) 


0.54 (0.56) 



7. General structure of program SIMUB 

The package for generation of production of i?-mesons and their decays, 
SIMUB, is kept under the directory SIMUB_V_X (X is the version number) 
which has the following main parts: 

• bb_gen - routines needed to generate bb events (FORTRAN, PYTHIA, 
HBOOK); 

• bb_f rg - routines performing string fragmentation and generation of 
5-mesons (FORTRAN, PYTHIA, HBOOK). The program is a part 
of BB_dec program but may be used as independent program; 
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• BB_dec - routines performing B-decays (C++, FORTRAN, PYTHIA, 
HBOOK, ROOT) and storing tlie results into standard HEPEVT- 
format Ntuple for further usage in the CMS detector simulation or 
in format JETSET for analysis; 

• include - a collection of common blocks for bb_gen, bb_frg and 
BB_dec. 

To simulate the S-meson production in pp coUisions with PYTHIA [2] 
(programs bb_gen and bb_f rg) we used the modified routines from the 
FORTRAN-based package described in [17]. 

The data flow between three main parts (bb_gen, bb_f rg, and BB_dec) 
of the package SIMUB_V_X is organized in the following two steps: 

• At the flrst step the 66-events are generated by the program bb_gen 
and written to Ntuple bb_X_YYY. 

• At the next step the program BB_dec reads 66-events from Ntuple 
bb_X_YYY and performs the string fragmentation (via call bb_f rg) and 
decay of B mesons. The decay modes and mechanisms are defined by 
user. The full information about B decay events is stored in Ntuple 
BB_dec.ntpl. Its format (HEPEVT or JETSET) is defined by user. 

8. Inheritance of classes in the program BB_dec 

The example of class inheritance of BB_dec program is shown in Fig. 2 
for the case of time-depended decay channel 

B\t),^{t) ^ v,i^ in-)V2{^ P+P-) , 

where Vi and V2 are the vector particles, l'^ and l~ are the leptons, P+ 
and P~ are the pseudoscalar mesons. The thick and thin arrows show 
the directions of class inheritance from mother to daughter classes. Thick 
arrows correspond to main inheritance ways. The name of each class co- 
incides with the name of subdirectory where this class is kept. The upper 
dashed-line box corresponds to the directory src. The lower dashed-line 
boxes unite the names of classes which are in the same subdirectory at the 
same level. In case when several classes belong to the same subdirectory, its 
name coincides with the name of the class where the thick incoming arrow 
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TEST-main T_Decay 



T_Utility T_Particle T_B_Event 



T_BPhysParain 



T TimeEvolution 



,T Baala2bb1b2 



T Baa1a2bb1b2Param 



T RandomBB 



T_CKM 



l 

T BOsI 



T Observ BVV 



T_BptT| 
T BOdI 



T GeneratorParam 



T_EvtLoopDecay , 



T_B_VV 
I 

T_V2_pri 
i 

T_B_VIIVp^ 

Figure 2: The class inheritance of the BB_dec program in case of simulation of 
decays B^{t),B^{t) Vi{-^ Z+rjFaH P^P~) 



shows to. The directories are included one into another in the directions 
opposite to thick arrows. 

The main directory src, shown in Fig. 2 as upper dashed-line box, con- 
tains five subdirectories T_Decay, T_Pcirticle, T_B_Event, T_Utility and 
TEST-main. The subdirectory T_Decay contains definition of the mother 
class with the same name as well as subdirectories with classes which are 
common for different decay channels. The subdirectory T_Particle con- 
tains the definition of the mother class with the same name and subdi- 
rectories with daughter classes T_BOs, T_Bpm and T_BOd which define the 
properties and decay modes for B*^ (-8°), B"^ and 5° (5°) mesons, respec- 
tively. The directory T_B_Event contains the definition of the mother class 
with the same name and daughter class T_EvtLoopBDec with a loop over 
all events with S-mesons. In the loop, each event is read from the ROOT 
Tree bb_frg.root and, then, the decays of S-mesons are performed ac- 
cording to modes defined in the samples of classes T_BOs, T_Bpm and T_BOd 
both for particles and anti-particles. The directory T_Utility consists of 
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the auxiliary classes and functions. The directory TEST-main contains the 
main programs for testing of classes. 

The directory src also includes the file T_EvtLoopBdec_main. C con- 
taining a mainO function which is not shown in Fig. 2. In this file, user 
can find the example of loop over the events. 

9. Monte Carlo methods for generation of decays 

B%t)^Vi{^ in')V2{^ P+P-) 

In class T_VertexB_VllVpp, the variables describing the decay channel 
and MC method for random generation arc defined. The member func- 
tion T_VertexB_VllVpp: :DecayB_VllVpp generates four random numbers 
- proper decay time t and angular variables cos©;,j, cosOai and x ~ a-c- 
cording to distribution function given by Eq. (8). Then, DecayB_VllVpp 
calculates Lorenz vector characterizing the decay vertex in Laboratory sys- 
tem. 

Two MC methods of random numbers generation are implemented in 
the class T_VertexB_VllVpp. 

The first method is based on the filling of the large single-dimensional 
array f 4_Integ [n4_cells] of real numbers which represents numerically 
the four-dimensional distribution function corresponding to /o(©ai ; ©fei ; X; t) 
given by Eqs. (8) and (11). The array f4_Integ contains n4_cells = 
f Npx X f Npy X fNpz X fNpt elements, where the variables f Npx, f Npy, 
fNpz and fNpt define the MC generator resolutions. The array f4_Integ 
is filled in the constructor of class T_VertexB_VllVpp where the mem- 
ory space for the array is reserved during all time of the existence of 
this class sample. The array f4_Integ is used in the member function 
T_VertexB_VllVpp : :GetRandom4 performing fast generation of the four 
random values of cos 6;,^, cos G^j, x t ^■ 

The second method involves sequential generation of random numbers 
according to approach which is based on usage of single-dimensional dis- 
tribution functions given by Eqs. (10) and (12). This method was already 
deseiilxxl in section 1. In this case the constructor T_VertexB_VllVpp does 

^In the function GetRandom4 we used the algorithm which is analogous to that was 
realized in the ROOT class TF3 [18] for generation of three random numbers distributed 
according to 3-dimensional probability function. 
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not fill any arrays and, therefore, it does not reserve the memory and may 
be used for definition of big number of particles decaying in the same mode. 

These methods have been tested in case of MC simulation of the "golden' 
decay mode Jli^{-^ K+K') with values of 

physics parameters fixed according to section 6. Both methods demon- 
strated the identical results with suitable time and memory consuming. 

The time-dependent angular correlations in the decays of type B^{t) 
Vi(— s> l^l^)V2{^ P^P~) can be used in data analysis to extract CP- 
violating asymmetries, — B mixing parameters and other observables 
[13, 19, 20]. Therefore, the decay dynamics described by angular correla- 
tions (8) should be included without fail into MC generators developed for 
physics studies. 




Figure 3: Comparison of time and angular distributions for decay 
^ J/^(^ n+n-)(t){-^ K+K-) generated by SIMUB and 
PYTHIA. For time t we use units Imm/c 3.33 x 10~^^ sec 
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In Fig. 3 wc compare time and angular distributions in the helicity 
frame for decay B'^{t), Bg{t) J/ip{^ ii'^ 4>{-^ K^K^) generated 
by SIMUB and PYTHIA. The difference between these two generators 
becomes more clear when slices over time and angular variables are cut 
from the kinematical phase space. The angular distributions shown in 
Fig. 3 are built for the following kinematical shces: 
X ■ 0.01 mm/c <t < 0.19 mm/c; 
cos©^+ : 0.01 mm/c <t< 0.19 mm/c, 0.3 rad < x < 0.6 rad; 
cos©x+ '■ 0.01 mm/c <t< 0.19 mm/c, 0.3 rad < x < 0.6 rad, 
0.7 < cos < 1 . 

In case of PYTHIA usage. Fig. 3 shows the uniform distributions for 
angular variables cos Qbi , cos Qai ^-nd x because of lack of time-dependent 
angular correlations. Due to this simple reason PYTHIA can not be used 
for MC studies of dynamics of sequential two body decays of B mesons 
in the channels with intermediate vector mesons. Unfortunately, due to 
various technical reasons, the both other well known packages, QQ [3] and 
EvtGcn [4], turn out to be also not suitable for study of angular correlations 
in decays of the type B — > Vi(/x"'"/x~) V2{P^P~). 

Conclusion 

General scheme of program chain and data flow for MC simulation of 
S-meson production and decays have been developed and realized in the 
package SIMUB which is now in progress. This scheme was tested for 
case of semileptonic decays B D^*^ jiv and sequential decays of the type 
B Vi{^ fi^/j^)V2{^ P^P^). The SIMUB package provides imique 
tools for MC physics studies of dynamics of semileptonic and "golden" B- 
decay channels with taking into account B^ — ^ oscillations and angular 
correlations. 

The starting version of the package SIMUB _V_0 is installed and tested at 
CERN on Linux platform. The package and documentation for user one 
can find on the SIMUB Home Web Page [1] . There is link to this page from 
CMS -B-Physics Croup Web Page: http://cmsdoc.cern.ch/~bphys. Now 
the SIMUB is used within CMS collaboration for exclusive S-trigger stud- 
ies. 
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